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This  paper  is  concerned  with  the  solution  of  a  real  world  Bilevel  Linear 
Program  (BLP)  of  the  highway  Network  Design  Problem  (NDP),  based  on  empirical 
data  for  Tunisia.  The  problem  includes  2,683  variables  (2,571  at  the  lower 
level)  and  820  constraints;  it  is  far  beyond  the  capacity  of  any  existing 
algorithm.  A  new  algorithm  is  devised  and  the  solution  is  found  successfully 
despite  the  NP-hardness  of  BLP  and  the  large  scale  of  the  problem.  The  computa- 
tional difficulties  of  BLP  should  not  present  a  barrier  against  the  effective 
use  of  such  a  powerful  optimization  technique. 


Bilevel  Linear  Programming  (BLP)  is  a  mathematical  model  of  an  organiza- 
tional hierarchy  in  which  two  decision  makers  have  to  choose  their  strategies 
from  the  polyhedron  {(x,y):  Ax  +  By  <  b;  x,  y  >  0).  The  upper  decision  maker, 
who  has  control  over  x,  makes  his  decision  first,  hence  fixing  x  before  the 
lower  decision  maker  selects  y  [Bialas  and  Karwan  1982  and  1984,  Bard  1983, 
and,  Ben-Ayed  1988] . 

The  ability  of  BLP  to  represent  decentralized  decision-making  problems 
motivated  the  development  of  several  applications  of  the  model,  including 
resource  allocation  [Cassidy  1971,  Fortuny-Amat  and  McCarl  1981],  military 
planning  [Bracken  and  McGill  1973,  1974a  and  1974b,  Shere  and  Wingate  1976, 
Falk  1977]  and  government  policy  [Candler  and  Norton  1977,  DeSilva  1978,  Bialas 
et  al.  1980,  Falk  and  McCormick  1982,  Candler  and  Townsley  1982  and  1985]. 
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Recently,  attention  has  been  focused  on  applying  bilevel  formulations  to  a 
transportation  decision  making  problem  referred  to  as  the  Network  Design 
Problem  (NDP)  [LeBlanc  and  Boyce  1986,  Marcotte  1986,  Ben-Ayed,  Boyce  and  Blair 
1988].  NDP  is  concerned  with  the  optimal  improvement  of  a  transportation 
network  in  order  to  minimize  the  system  total  costs. 

Most  BLP  applications  are  limited  to  the  theoretical  formulation,  thereby 
lacking  the  practical  significance  of  solving  real  problems.  Actually,  the 
overwhelming  majority  of  real-world  problems  are  formulated  and  solved  as 
single-level  (single-objective)  programs  even  when  they  are  virtually  bilevel. 
This  is  mainly  due  to  the  inefficiency  associated  with  BLP  algorithms  [Ben-Ayed 
and  Blair  1988],  which  presents  a  real  barrier  against  the  effective  use  of  the 
model.  In  trying  to  motivate  further  applications  of  BLP  going  beyond  the 
theoretical  formulations  and  the  small  illustrative  examples  to  interesting 
real  world  problems,  Ben-Ayed,  Blair  and  Boyce  [1988]  presented  a  BLP  formula- 
tion of  the  interregional  highway  NDP  based  on  empirical  data  for  Tunisia.  This 
paper  gives  the  solution  to  the  problem  formulated  there.  This  problem  includes 
2,683  variables  (2,571  of  which  are  lower)  and  820  constraints,  which  makes  it 
too  large  to  be  solved  by  any  of  the  existing  algorithms. 

Section  1  presents  the  problem.  In  Section  2,  we  develop  a  new  algorithm 
that  takes  advantage  of  the  special  structure  of  the  problem.  There  are  two 
reasons  for  having  a  BLP  formulation;  first  the  user-optimized  flow  requirement 
(user-equilibrium),  and  second  the  nonconvex  improvement  functions.  The 
algorithm  deals  with  each  of  the  two  lower  problems  separately;  at  each 
iteration  we  try  to  find  a  better  compromise  with  the  user,  while  including  the 
smallest  possible  number  of  nonconvex  improvement  functions  to  get  the  exact 
solution  with  the  minimum  computation  effort.  In  Section  3,  the  solution  is 
presented  and  the  results  are  analyzed.  We  conclude  in  the  last  section  with 
some  suggestions. 


1.  The  Problem  to  be  Solved 

In  this  section  we  give  a  short  description  of  the  formulation  resulting 
from  the  empirical  study  conducted  by  Ben-Ayed,  Blair  and  Boyce  [1988]  on  the 
Tunisian  interregional  highway  transportation  network.  We  are  concerned  with  a 
system  of  roads  connecting  the  regions  of  Tunisia.  Each  region  is  represented 
by  a  centroid-node  which  is  usually  the  largest  city  in  the  region.  The  traffic 
entering  and  leaving  a  centroid-node  is  assumed  to  be  that  of  the  whole  region; 
every  centroid-node  is  simultaneously  an  origin-node  and  a  destination-node. 
The  other  nodes  in  the  network  are  intersections  of  roads  generally  correspond- 
ing to  intermediate  cities.  A  road  connecting  any  two  cities  is  called  a  link 
and  a  sequence  of  links  is  called  a  route.  We  have  empirical  data  on  the  number 
of  people  who  want  to  travel  from  each  region  to  each  other  region.  Typically, 
a  traveller  has  a  choice  of  routes  and  makes  his  decision  based  on  road 
conditions.  Improvements  in  various  links  in  the  network  sill  cause  some 
travellers  to  change  their  choice  of  routes.  We  will  not  provide  details  about 
the  empirical  data,  the  construction  of  the  objective  functions  and  the 
constraints  in  the  formulation;  however,  a  brief  discussion  on  the  inclusion  of 
the  bilevel  presentation  is  relevant  to  the  object  of  this  paper. 

The  principal  variables  in  the  formulation  are  the  traffic  flow  and  the 
added  capacity  resulting  from  the  improvement,  both  measured  in  Passenger  Car 
Units  per  hour  (PCU).  The  designers  of  a  road  network  wish  to  find  the  best 
combination  of  improvements  to  make.  Their  objective  takes  into  account  the 
cost  of  improvements  and  also  factors  which  depend  on  the  pattern  of  road  use, 
including  system  travel  time,  system  operating  cost  and  frequency  of  accidents. 
The  users  of  the  network  decide  which  routes  to  use  primarily  based  on  their 
individual  cost;  therefore,  their  choices  do  not  necessarily  coincide  with  the 
choices  that  are  optimal  for  the  system  [Wardrop  1952].  The  system  can  in- 
fluence users   choices,  however,  by  improving  some  links  and  making  them  more 
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attractive  than  others.  When  deciding  on  these  improvements,  the  system  tries 
to  influence  the  users'  preferences  in  order  to  minimize  total  system  costs. 
Thus  in  our  BLP,  the  upper  decision  maker  chooses  the  amount  of  improvement  to 
make  on  each  of  the  different  roads  (links  in  the  network)  while  the  lower 
decision  maker  chooses  how  much  traffic  there  will  be  on  each  road.  A  similar 
approach  is  used  by  LeBlanc  and  Boyce  [1986],  and  Ben-Ayed,  Boyce  and  Blair 
[1988] . 

In  addition  to  the  user-optimized  flow  requirement,  BLP  formulation  allows 
the  analysis  of  nonconvex  piecewise  linear  functions  [Ben-Ayed,  Blair  and  Boyce 
1988].  A  typical  improvement  function  might  be: 

f(Z)  =  10Z  0  <  Z  <  200 

=  2000  +  3(Z-200),   200  <  Z  <  500  (1) 

=  2900  +  8(Z-500),   500  >  500. 

where  Z  is  the  amount  of  improvement  of  a  link  and  f(Z)  is  the  cost  of  that 
improvement.  Note  that  f(Z)  is  neither  convex  nor  concave,  thus  making  it 
difficult  to  analyze.  However,  the  function  can  be  included  in  the  upper 
objective  of  a  BLP  with  an  appropriate  lower  objective.  The  technique  is  based 
on  the  following: 

f(Z)  =  F(Z,W)  =  10Z  -  7WX  +  5W2,  where  W^MAXfO,  Z-200}  ,  W2=MAX(  0 ,  Z-500]  . 

The  variables  Wm  are  controlled  by  the  lower  problem;  the  lower  objective 
trying  to  make  the  Wm  as  small  as  possible,  and  the  constraints  W^  >  Z-200,  W2 
>  Z-500  guarantee  that  Wm  are  least  upper  bounds.  This  technique  can  be  used 
for  any  piecewise  linear  function. 

Two  drawbacks  should  be  mentioned.  First,  the  method  described  does  increase 
the  number  of  lower  variables,  and  this  increase  may  be  considerable  if  there 
are  many  piecewise  linear  functions  in  the  problem.  Second,  the  resulting  BLP 
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is  typically  distinguished  by  the  conflict  between    its  upper  and  lower 
objectives;  there  are  many  variables  which  the  upper  objective  wants  to  be  as 
large  as  possible,  while  the  lower  objective  wants  them  as  small  as  possible. 
Experience  suggests  this  type  of  BLP  is  difficult  to  solve. 

Including  user-optimized  flow  and  nonconvex  improvement  functions,  the  BLP 
problem  obtained  by  Ben-Ayed,  Blair  and  Boyce  [1988]  is  synthesized  in  the 
following: 

MINZ  F(Z,W)  +  G(C) 

where  [C,  C,  W,  X  and  Y]  solve: 

MIN  F(W)  +  G(C) 

subject  to: 

%(¥)  =  0 

(2) 

H2(X,Y)  =  0 

H3(Z,C,X)  >  0 

H4(Z,C,X)  >  0 

H5(Z,W)  >  0 

H6(Z)  >  0 

Z,  C,  C,  W,  X,  Y  >  0 

where  the  function  F  is  the  improvement  cost.  The  function  G  is  the  inter- 
regional system  travel  cost.  The  vector  X  is  the  amount  of  traffic  on  each 
link.  Any  traveller  from  one  origin-node  to  a  destination-node  has  a  choice  of 
different  routes  to  use,  corresponding  to  different  sets  of  links.  We  use  the 
components  of  the  vector  Y  to  indicate  the  number  of  people  using  each  route. 
These  numbers  determine  the  numbers  of  users  of  each  link,  shown  by  the  vector 
X;  the  relationship  between  X  and  Y  is  given  by  the  equation  constraints 
H2(X,Y)  while  the  requirement  that  specified  numbers  of  people  go  between 
specified  nodes  is  given  by  the  constraints  H^Y)  referred  to  as  the  conserva- 
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tion  of  flow  conditions.  F  is  the  lower  objective  of  nonconvex  improvement  cost 
functions.  The  components  of  the  vector  Z  are  the  amounts  of  improvement  on 
each  of  the  different  links.  W  is  the  vector  of  the  second-stage  variables  used 
in  computing  the  improvement  cost  functions  as  described  previously  in  (1).  G 
gives  the  interregional  travel  costs  on  each  link  for  the  user.  C  and  C  are 
vectors  of  second-stage  variables  giving  total  travel  costs  on  each  link  for 
the  system  and  the  user,  respectively.  H3(Z,C,X)  and  H^(Z,Z,X)  are  the  ine- 
quality constraints  needed  to  the  formulation  of  convex  travel  cost  functions 
and  the  effect  on  improvement  on  decreasing  those  costs  for  the  system  and  the 
user,  respectively.  Hr(Z,W)  correspond  to  the  nonconvex  piecewise  linear 
functions  constraints.  Hg(Z)  are  the  improvement  limit  constraints.  A  more 
detailed  description  of  the  variables  and  the  functions  is  provided  in  the 
Appendix  A. 

In  the  above  formulation  (2),  the  total  number  of  variables  (see  Appendix  A) 
is  2,683  and  the  total  number  of  constraints  is  820,  excluding  upper  and  lower 
bound  constraints  (improvement  limit  constraints  and  nonnegativity  const- 
raints). Those  numbers  are  partitioned  as  follows: 

Variables  Constraints 

(excluding  bound  constraints) 

112  (X) 

112  (C)  120  (conservation  of  flow  conditions  ) 

112  (C)  112  (link  flow  definition  constraints) 

112  (Z)  224  (system  travel  costs  formulation) 

104  (W^)  224  (cumulative  user  costs  formulation) 

36  (Wo)  68  (convex  improvement  costs  formulation) 

2,095  (Y)  72  (nonconvex  functions  formulation) 


2,683  820. 

The  Z  are  the  only  variables  controlled  by  the  upper  problem;  therefore 
among  the  2,683  variables,  2,571  are  lower.  The  nonconvex  improvement  cost 
functions  all  have  3  pieces;  two  constraints  are  needed  by  each  link.  The 
convex  improvement  functions  require  only  one  constraint  per  link  because  they 


all  consist  of  two  pieces  only. 

An  ideal  solution  (lower  bound)  for  any  BLP  minimization  problem  can  be 
obtained  by  ignoring  the  lower  objective  function  and  solving  the  problem  as  an 
LP.  If  we  call  (Z°,  X°,  C°,  C°,  Y°,  W°)  the  obtained  solution,  the  ideal 
solution  is  equal  to  F(Z  ,W  )  +  G(C  ).  To  get  an  incumbent  solution  (upper 
bound),  we  substitute  Z  by  Z  and  solve  the  lower  LP  problem.  If  we  call  (X  , 
C  ,  C  ,  Y  ,  W  )  the  obtained  solution,  the  incumbent  solution  is  equal  to 
F(Z  ,W  )  +  G(C  ).  Let  us  refer  to  the  gap  between  incumbent  and  ideal  solutions 
as: 

gap  =  (incumbent  -  ideal)/ideal . 

The  gap  for  the  BLP  problem  (2)  is  larger  than  111  %  [Ben-Ayed  1988];  this  high 
value,  mainly  due  to  the  nonconvex  functions  formulation,  makes  the  task  of 
heuristic  algorithms  extremely  difficult,  if  not  impossible.  Heuristic  al- 
gorithms, such  as  the  Parametric  Complementary  Pivot  algorithm  (PCP)  [Bialas- 
Karwan-Shaw  1980,  and  Bialas-Karwan  1984]  and  the  Grid  Search  Algorithm  (GSA) 
[Bard  1983],  lead  to  very  unsatisfactory  results  when  applied  to  BLP  problems 
with  combinatorial  nature,  such  as  Knapsack  problems  or  nonconvex  formulation 
problems  [Ben-Ayed  and  Blair  1988]. 

The  problem  is  even  harder  for  branch  and  bound  or  related  implicit  enumera- 
tion techniques  [Falk  1973,  Gallo  and  Ulkiicu  1977,  Fortuny-Amat  and  McCarl 

1981,  Bialas  and  Karwan  1982,  Papavassilopoulos  1982,  Candler  and  Townsley 

1982,  Bard  and  Falk  1982,  Bard  and  Moore  1987].  The  Bard-Moore  algorithm,  which 
is  supposed  to  be  one  of  the  most  efficient,  does  not  accept  a  problem  with 
more  than  100  lower  variables;  our  problem  has  more  than  25  times  this  number. 

As  for  any  large  NP-hard  problem,  the  use  of  general  algorithms  is  not 
efficient.  A  new  algorithm  taking  advantage  of  the  special  structure  of  the 
problem  has  to  be  devised  to  solve  the  problem  at  hand. 


(4) 


2.  A  Proposed  Algorithm 

Formulation  (2)  shows  that  there  are  two  separate  lower  objective  functions. 
In  fact,  each  of  the  two  objectives  has  its  own  constraints  that  are  indepen- 
dent of  the  constraints  of  the  other  objective,  provided  that  Z  are  fixed  (they 
are  fixed  by  the  upper  objective).  This  can  be  better  seen  when  the  BLP  problem 
(2)  is  restated  as  follows: 

MINZ>0  F(Z,W)  +  G(C)  (3) 

such  that  {W)  solve  (4): 

MIN  F(W) 

subject  to: 

H5(Z,W)  >  0 

W  >  0 

and  {C,  C,  X  and  Y)  solve  (5): 

MIN  G(C) 

subject  to: 

H(X,C,C,X,Y)  >  0 

C,  C,  X,  Y  >  0 

where  H(X,C,C,X,Y)  >  0  substitutes  the  constraints  H1(Y)  =  0,  H2(X,Y)  =  0, 
H3(Z,C,X)  >  0,  H4(Z,C,X)  >  0,  and  H6(Z)  >  0. 

The  first  lower  problem  (4)  corresponds  to  the  nonconvex  functions  formula- 
tion and  the  other  lower  problem  (5)  corresponds  to  the  equilibrium  problem. 
Because  of  the  NP-hardness  of  BLP,  it  is  more  convenient  to  take  advantage  of 
the  structure  of  the  problem  and  solve  it  as  two  separate  smaller  problems. 

As  we  mentioned  in  our  discussion  of  improvement  cost  functions  in  the 
previous  section,  the  BLP  corresponding  to  (4)  is  difficult.  However,  the 
situation  is  different  for  the  BLP  corresponding  to  the  equilibrium  problem 
(5);  the  travel  cost  functions  of  the  system  and  the  user  are  very  similar 


(5) 
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[Ben-Ayed,  Blair  and  Boyce  1988].  The  coefficients  of  the  two  functions  always 
have  the  same  sign,  and  are  even  fairly  proportional;  therefore  the  use  of  BLP 
heuristics  is  safe  for  this  second  part. 

Our  approach  to  satisfy  (5)  tries  to  identify  a  solution  that  is  good  for 
both  the  upper  and  lower  objectives  by  using  a  convex  combination  of  the  two: 

tG(C)  +  (l-t)G(C) 

so  that  we  are  looking  at  the  problem: 

MINZ  F(Z,W)  +  tG(C)  +  (l-x)G(C) 

such  that  {Z,  C,  C,  X  and  Y)  solve: 

H(X,C,C,X,Y)  >  0 

Z,  C,  C,  X,  Y  >  0 

and  such  that  (W)  solve:  (6) 

MIN  F(W) 

subject  to: 

H5(Z,W)  >  0 

W  >  0. 

In  general,  a  BLP  solution  (x  ,y  )  is  called  feasible  if  it  belongs  to  the 
polyhedron  (Ax°  +  By  <  b;  x  ,  y°  >  0};  and  satisfies  the  optimality  of  the 
lower  solution  y  to  the  lower  problem  when  the  upper  variable  x  is  fixed  at 
x  .  In  our  case,  for  a  feasible  solution  to  (6)  to  be  feasible  to  (2),  it  must 
be  the  case  that,  for  the  given  choice  of  the  upper  variables  Z,  the  values  of 
the  lower  variables  are  such  as  to  minimize  the  lower  objective  G.  This  will 
clearly  be  the  case  for  t  =  0.  To  make  the  upper  objective  as  small  as  pos- 
sible, we  want  to  find  the  largest  x  for  which  this  is  the  case. 

We  will  assume  that  if,  for  some  T  ,  the  solution  to  (6)  is  feasible,  then 
the  same  will  be  true  for  any  0  <  t  <  t  .  This  is  not    always  true,  but  counter- 
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examples  are  not  easy  to  construct,  and  a  simulation  suggests  this  is  a  safe 
assumption.  Granted  this  assumption,  the  optimal  t,  call  it  x    ,  can  be  located 
by  a  binary  search,  which  is  much  quicker  than  a  grid  search,  as  used  in  Bard 
[1983]. 

For  a  given  t,  (6)  is  still  a  BLP  with  too  many  lower  variables  to  be  solved 
in  a  straightforward  way.  We  approach  this  problem  by  trying  to  identify  links 
on  which  no  improvement  will  be  made.  For  those  links,  we  replace  the  nonlinear 
improvement  cost  function  f(Z)  by  aZ  +  B,  chosen  so  that  aZ  +  /3  <  f(Z)  for  all 
Z.  We  are  temporarily  assuming  that  the  costs  for  those  links  we  think  are 
unlikely  to  be  improved  are  less  than  they  really  are.  This  produces  a  smaller, 
more  tractable,  problem.  If  the  solution  to  this  smaller  problem  does  not  use  a 
link  with  the  simplified,  cheaper  cost  it  certainly  would  not  want  an  improve- 
ment with  the  real  cost.  In  the  case  where  the  small  problem  does  recommend  an 
improvement  based  on  the  simplified  cost,  we  must  solve  a  new  problem  with  the 
real  cost  as  part  of  the  formulation.  Thus  we  solve  a  sequence  of  problems  in 
which  links  with  real  costs  are  added  when  it  seems  necessary,  but  not  other- 
wise. Fortunately,  in  practice  the  process  terminates  (no  new  links  added)  with 
a  relatively  small  number  of  links  showing  improvement. 

Let  us  call  S  the  set  of  the  links  having  convex  or  linear  cost  improvement 
cost  functions,  T  the  set  of  the  links  having  nonconvex  improvement  cost 
functions,  T  the  subset  of  T  containing  all  links  for  which  Z  are  different 
from  0,  and  T  the  subset  containing  the  remaining  elements  of  T.  The  links 
belonging  to  T  are  considered  as  candidates  for  improvement  and  must  have 
their  three-piece  improvement  cost  functions  f(Z)  included  in  the  formulation. 
Let  (Z1,  X1,  C1,  C1,  Y1,  W1)  be  the  solution  of  the  LP  where  all  nonlinear 
functions  f(Z)  in  (6)  are  substituted  by  aZ  +  J3 ,  for  each  link  of  T.  We  obtain 
the  following  BLP  problem  with  relatively  few  lower  variables: 
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MINZ  ZagS  Fa(Za,Wma)  +  Za£T1  Fa(Za,Wma)  +  Za£T2  (aQZa+l3a)  +  xG(C)  +  (l-x)G(C) 
such  that  [Z,  C,  C,  X  and  Y}  and  {Wma,  m=l,Ma,  aeS}  solve: 

H(X,C,C,X,Y)  >  0 
Z,  C,  C,  X,  Y  >  0 

H5a(Za>Wma)  *  ° 

Wma  *  °  (7) 

and  such  that  {Wma,  m=l,Ma,  acT1}  solve: 

MIN  EaeTl  Zm=l,Ma  Wma 
subject  to: 

Wma  "  Za  -  "^ma 
Wma  *  °- 

where  Ma  is  the  number  of  breakpoints  in  the  piecewise  improvement  cost 
function  of  link  a;  Fa  is  the  upper  objective  function  of  nonconvex  improvement 
cost  on  link  a;  and,  Hca  gives  the  corresponding  constraints. 

Let  us  call  (Z2,  X2,  C2,  C2,  Y2,  W2)  the  solution  of  (7).  If  there  exists  a 

1        ") 
link  j  such  that  j  does  not  belong  to  T   and  Z-   is  different  from  0,  then  j 

has  to  be  subtracted  from  T  ,  added  to  T  ,  and  (7)  has  to  be  resolved  until  all 

?  1 

links  a  having  Za   positive  are  elements  of  T  ,  in  which  case  the  exact 

solution  of  (6)  is  found. 

We  have  constructed  from  the  main  BLP  problem  (2)  two  other  secondary  BLP 
formulations:  a  parameterized  BLP  (6)  and  a  BLP  with  few  lower  variables  (7). 
Problem  (7)  is  solved  to  give  a  feasible  and  optimal  solution  to  (6).  If  this 
solution  is  user-optimized  (i.e.  solves  the  lower  LP  problem  (5)  when  Z  is 
fixed  at  its  optimal  value  for  (7)),  then  it  is  feasible  for  (2)  and  a  higher 
value  of  t  can  be  tried.  Otherwise,  we  continue  our  search  for  the  optimum  by 
decreasing  the  value  of  t. 

For  any  solution  (Z1,  X1,  C  ,  C  ,  Y1,  W1),  we  have  to  distinguish  between 
two  different  values.  The  first  value  is  that  of  the  upper  objective  function 
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of  formulation  (2),  equal  to  F(Z,W)  +  G(C);  let  us  call  it  01.  The  second  value 
is  the  actual  cost  incurred  by  the  system  given  a  specific  value  Z1  of  Z;  let 
us  call  it  A1.  A1  includes  actual  improvement  cost  and  actual  system  travel 
cost.  The  actual  improvement  cost  is  obtained  by  plugging  Z1  in  the  actual 
improvement  functions  f(Z),  whether  they  are  linear,  convex,  or  nonconvex.  The 
actual  travel  cost,  on  the  other  side,  is  obtained  after  we  solve  the  lower  LP 
problem  (5)  with  Z  fixed  at  Z1,  obtain  the  solution  (XJ,  C  J ,  C  J ,  Y  J )  and  plug 
the  XJ  into  the  convex  system  travel  costs;  in  other  words,  A1  is  equal  to 
fCZ1)  +  G(CJ).  The  first  value  01  corresponds  to  a  potential  optimal,  but  not 
necessarily  feasible  solution;  it  allows  to  find  ideal  solutions  to  the  BLP 
problem  (2).  The  second  value  A1  corresponds  to  a  feasible,  but  not  necessarily 
optimal,  solution;  it  allows  to  find  incumbent  solutions  for  (2). 

The  proposed  algorithm  is  an  iterative  procedure  that  tries  at  each  itera- 
tion to  reduce  the  gap  between  ideal  solution  and  incumbent  solution.  The 
procedure  terminates  when  the  gap  is  brought  below  a  desired  value  e,  or  when 
the  number  of  iterations  exceeds  a  fixed  limit  N.  The  following  flow  chart 
gives  a  general  description  of  the  iterative  procedure;  more  details  are 
provided  by  Appendix  B: 

Searching  for  x   (Using  Binary  Search) 


Start  with  t=1 

- 

A.  Solve  (6) 

B.  Test  whether  solution  of  (6) 
is  feasible  to  (2) 
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The  most  difficult  step  in  the  algorithm  is  step  A,  in  which  (6)  has  to  be 
solved.  Using  the  concepts  explained  above  in  this  section,  the  solution  of  (6) 
is  obtained  as  shown  in  the  following  flow  chart: 

t  Fixed 


Al.  Start  with  TX=0,  T2=T 

A2.  Solve  (7) 

A3.  Test  whether  new  elements 
should  be  added  to  T 

New 

elements  added 

No  element? 

;  added 

A/.      M~r»   Tl     r^r.T   T2 

K      (  C\        r.~1,r~^ 

J.   ,   new 

Several  improvements  can  be  added  to  the  algorithm  to  make  it  more  effi- 
cient. The  starting  ideal  and  incumbent  solutions  can  be  found  in  a  more  effi- 
cient way.  A  better  ideal  solution  for  (2)  can  be  obtained  when  solving  the  LP 
obtained  after  substituting  at  x  =  l  all  nonlinear  functions  f(Z)  in  (6)  by  ccZ+B; 
there  is  no  situation  better  for  the  system  than  having  cheap  improvement  and 
system-optimized  flow.  Also,  the  actual  cost  incurred  by  the  system  when  the 
problem  is  solved  as  a  BLP  must  be  lower  than  the  actual  cost  incurred  when  a 
simpler  technique  is  used;  otherwise,  there  is  no  justification  for  using 
complex  techniques.  Therefore,  to  obtain  an  incumbent  solution,  we  can  just 
solve  the  problem  as  an  LP;  we  approximate  improvement  cost  functions  by  their 
best  linear  fit,  ignore  the  lower  objective  and  solve  the  problem.  Let  us  call 
(Z°,  X°,  C°,  C°,  Y°)  the  obtained  solution,  A0  is  the  value  of  the  starting 
incumbent  solution.  Using  those  techniques  makes  the  gap  between  starting  ideal 
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solution  and  starting  incumbent  solution  as  low  as  3  %  instead  of  the  111  % 
found  earlier  in  the  previous  section. 

Moreover,  since  the  procedure  is  iterative,  it  is  more  efficient  to  use  the 
information  from  the  previous  iterations  in  order  to  make  the  next  one  easier 
to  solve.  One  way  to  do  this  is  to  introduce  to  the  problem  the  constraint: 

L  <  F(Z,W)  +  G(C)  <  U 

where  L  and  U  are  the  values  of  the  current  ideal  and  incumbent  solutions, 
respectively. 

Similarly,  the  solution  of  step  2  may  require  solving  (7)  more  than  once. 

9     9     9     9     9 

Let  us  assume  that  we  solved  (7),  obtained  the  solution  (Z  ,  X  ,  C  ,  C  ,  Y  , 

9 

W  )  and  have  to  add  a  new  element  and  resolve  (7).  We  use  the  fact  that  the  new 

cost  will  be  at  most  as  low  as  the  value  we  got  for  the  objective  when  we  used 

9  1 

cheap  linear  improvement  costs  (a-Z-   +  J3-)  for  the  links  j  to  be  added  to  T  ; 

which  gives  a  lower  bound  for  the  solution.  An  upper  bound  can  also  be  found  by 

9 

noticing  that  the  new  cost  is  at  least  as  high  as  the  actual  cost  f(Z  ) 

9 

incurred  by  the  system  for  an  added  capacity  equal  to  Z  .  The  difference 

9 

between  upper  and  lower  bounds,  when  step  2  has  to  be  repeated,  is  E^ti  (ai^i 
+  Bj)  -  f(Z2). 

Finally,  a  feasible  solution  to  a  BLP  problem  is  obtained  by  finding  the 
optimal  solution  to  the  lower  problem  after  substituting  values  for  the  upper 
variables.  In  our  case,  any  optimal  solution  (Z1,  X1,  C1,  C  ,  Y1,  W1)  of  the 
BLP  problem  (6)  satisfies  the  optimality  requirement  of  W1  to  the  lower  problem 
(A);  therefore  (Z1,  X1,  C1,  C1,  Y1,  W1)  is  a  feasible  solution  of  the  BLP 
problem  (2)  if  it  satisfies  the  optimality  requirement  for  the  second  problem 
(5).  Even  when  the  solution  (Z1,  X1,  C1,  C1,  Y  ,  W1)  is  not  feasible  a  new 
solution  (Zk,  Xk,  Ck,  Ck,  Yk,  Wk)  can  be  obtained,  with  feasibility  guaranteed. 
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It  is  known  from  the  theory  of  Linear  Programming  [see  for  instance  Chvatal 
1983,  Chapter  10  pp.  148-168]  that  given  two  LP  problems  LP1:  (MIN  ex,  Ax  =  b, 
x  >  0}  and  LP2:  (MIN  ex,  Ax  =  b' ,  x  >  0),  where  b  and  b'  are  different,  any 
basic  feasible  solution  x  of  LP2  that  has  the  same  nonbasic  variables  as  the 
optimal  solution  x  of  LP1,  is  an  optimal  solution  of  LP2.  Therefore  to 
generate  feasible  solutions  for  BLP  problem  (2),.  which  is  equivalent  to 
generating  optimal  solutions  to  the  lower  problem,  we  can  solve  BLP  problem  (7) 
with  the  additional  requirement  that  all  nonbasic  variables  of  the  solution 
(XJ,  CJ ,  C  J ,  Y  J )  of  the  lower  problem  (5)  with  Z  fixed  at  Z1  and  all  nonbasic 
slack  and  surplus  variables  be  locked.  This  is  added  to  the  algorithm  as  shown 
in  Appendix  B. 

Regarding  the  reliability  of  the  algorithm,  there  are  no  heuristics  involved 
in  finding  the  incumbent  and  the  ideal  solutions.  The  values  given  by  the 
algorithm  for  those  solutions  are  exact.  As  long  as  the  algorithm  converges, 
the  solution  it  gives  is  the  actual  optimum  within  the  accuracy  desired. 
However,  it  is  possible  that  one  finds  bad  examples  for  this  algorithm  where 
the  convergence  is  not  satisfactory;  in  such  a  case  the  algorithm  ends  with  a 
large  gap  between  upper  and  lower  bounds,  but  does  not  give  a  wrong  solution. 

3.  Solution  of  the  Real  World  Problem 

When  trying  to  solve  the  BLP  problem  (7)  we  faced  a  major  constraint,  since 
we  could  find  no  reliable  BLP  software  that  could  handle  problems  this  large. 
Fortunately,  (7)  can  be  transformed  into  an  equivalent  Mixed  Integer  Program- 
ming problem  (MIP),  with  the  number  of  zero-one  variables  in  the  latter  equal 
to  the  number  of  lower  variables  in  (7).  It  should  be  emphasized  that  the  use 
of  MIP  does  not  mean  that  it  is  easier  to  solve.  MIP  and  BLP  are  both  NP-hard 
problems  and  are  equally  difficult  to  solve;  however,  MIP  is  older  and  benefits 


16 
from  the  availability  of  good  software  to  solve  it.  One  such  package  is  the 
Zero/One  Optimization  Method  (ZOOM),  which  is  an  extension  of  the  XMP  Linear 
Programming  system.  ZOOM  and  XMP  are  both  written  in  Fortran  77  and  revised  in 
1986  by  R.  E.  Marsten  from  the  Department  of  Management  and  Information  Systems 
at  the  University  of  Arizona.  A  description  of  the  XMP  library  is  given  in 
Marsten  [1981].  The  subroutines  of  ZOOM  and  XMP  are  used  by  our  program 
whenever  an  LP  or  an  MIP  needs  to  be  solved.  The  computational  study,  to  be 
presented  in  this  section,  was  conducted  on  the  CRAY  X-MP/24  supercomputer  of 
the  National  Center  for  Supercomputing  Applications  at  the  University  of 
Illinois  at  Urbana-Champaign. 

First,  we  show  how  to  transform  BLP  (7)  into  MIP.  The  constraint: 

Wm  =  MAX{(Z  -  qm),  0) 

is  equivalent  to: 

either:  Wm  =  Z  -  qm,  or:   Wm  =  0 

which  can  be  equivalently  restated,  after  introducing  zero-one  variables  Vma, 
as : 

Wm  -  Z  <  -qm  +  MVm,   -Wm  +  Z  <  qm  +  MVm 

(8) 

Wm  *  "(l-V,   "Wm  *  ^-V 

where  M  is  a  large  number  (it  is  sufficient  to  make  M  larger  than  the  maximum 
possible  road  improvement). 

The  formulation  of  (7)  as  MIP  is  obtained  by  replacing  the  lower  problem  in 
(7)  by  the  constraints  (8)  for  each  link  a.  The  constraints  Wm  <  M(l-Vm)  are 
trivial  because  of  the  nonnegativity  of  Wm;  we  take  them  out,  and  obtain  the 
following  problem  to  be  solved  whenever  we  need  to  solve  (7): 
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«INZ  ZaeS  Fa(Za'Wma)  +  EaeTl  Fa^a'Wma)  +  ZaeT2  (°aZa+6a)  +  *G^   +  (l-t)G(C) 
such  that  (Z,  C,  C,  X  and  Y}  and  {Wma,  m=l,Ma,  aeS)  solve: 

H(X,C,C,X,Y)  >  0 

Z,  C,  C,  X,  Y  >  0 

H5a(Za>Wma)  *  0 


and  such  that  {Wma,  m=l,Ma,  aeT  )  solve: 


(9) 


"Wma  +  Za  "  MVma  *  ^ma 

Wma  "  Za  "  MVma  *  "^ma 

*ma  +  MVma  *  M 
Vma  e  (0,  1). 

The  technique  used  to  transform  (7)  into  MIP  is  different  from  the  technique 
used  by  Fortuny-Amat  and  McCarl  [1981].  The  use  of  this  latter  would  have 
doubled  the  number  of  zero-one  variables  (requiring  a  zero-one  variable  for 
each  lower  variable  and  a  zero-one  variable  for  each  lower  constraint). 

The  large  size  of  the  problem  (2)  is  mainly  due  to  the  high  number  (2,095) 
of  routes  included  in  the  formulation.  Since  not  all  those  routes  are  going  to 
be  used,  the  problem  can  be  simplified  by  keeping  the  most  attractive  ones  and 
discarding  the  rest.  We  will  solve  the  problem  with  the  reduced  number  of 
routes  (we  refer  to  this  problem  as  the  small  version  of  the  problem),  and  we 
solve  also  the  original  problem  with  all  routes  included;  next  we  will  compare 
the  solutions  of  the  two  problems.  The  attractive  routes  are  found  by  solving 
both  the  user  equilibrium  problem  (5)  and  the  system  equilibrium  problem  (same 
as  (5)  except  that  G(C)  is  replaced  by  G(C))  as  we  increase  the  values  of  the 
entries  of  the  trip  matrix  to  make  the  users  and  the  system  choose  their  next 
desirable  routes.  When  solving  the  system  equilibrium,  we  vary  also  the 
improvement  cost  because  of  its  effect  on  the  routes  choice.  The  accumulated 
number  of  routes  chosen  by  both  the  system  and  the  user  reaches  389  as  the  flow 
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fluctuates  between  the  original  trip  matrix  and  three  times  that  matrix.  Those 
389  routes,  averaging  more  than  three  routes  per  origin  destination  pair  are 
the  ones  to  be  included  in  the  small  version  of  the  problem  as  compared  to  the 
2,095  included  in  the  large  version. 

About  15.25  minutes  CPU  time  and  1.4  million  words  (64  bits  per  word) 
computing  storage  space  were  required  to  obtain  the  solution  of  the  small 
version  of  the  problem.  At  iteration  0,  we  compute  starting  ideal  and  incumbent 
solutions  using  the  ordinary  techniques  of  BLP;  the  obtained  gap  is  111.36  %; 
next,  we  use  the  techniques  proposed  in  Section  2  to  get  better  starting 
solutions.  At  the  same  iteration,  we  solve  the  problem  as  an  LP  and  we  use  the 
actual  cost  incurred  by  the  system  (namely  257,037.230)  as  a  new  incumbent 
solution.  Then  we  start  the  iteration  1  with  T=l,  TflAX13-'-  anc*  TMIN=^'  When  we 
find,  in  the  first  step  (A1-A2  in  the  flow  chart)  of  iteration  1,  the  list  of 
links  that  are  candidates  for  improvement,  we  use  the  solution  as  a  new  better 
ideal  solution;  the  gap  decreases  to  3.03  %.  In  the  second  step  (A3-A4-A2  in 
the  flow  chart),  the  MIP  (9)  is  solved  with  only  six  integer  variables.  The 
solution  obtained  shows  that  no  more  integer  variables  need  to  be  added;  a 
slightly  improved  ideal  solution  is  found  for  the  problem  (new  gap  equals  2.93 
%) .  The  feasibility  of  the  solution  is  tested  in  step  3  (B  in  the  flow  chart); 
the  flow  is  not  user-optimized.  We  generate  a  feasible  solution  using  the 
solution  obtained  when  testing  feasibility;  a  new  incumbent  solution  is 
obtained  and  the  gap  is  reduced  to  2.57  %  at  the  end  of  the  first  iteration.  At 
iterations  2  and  3  there  is  no  success  in  improving  the  incumbent  solution. 
Iteration  4,  where  x  =  .125,  i^x  =  -250  and  t^IN  =  0,  gives  another  slight 
improvement  at  the  end  of  step  3,  when  generating  a  new  feasible  solution.  The 
procedure  is  stopped  at  Iteration  5  because  the  solution  seems  to  be  reaching  a 
steady  state;  Iteration  5  gives  exactly  the  same  solutions  as  Iteration  4  in 
all  steps.  The  final  solution  retained  incurs  a  cost  to  the  system  equal  to 
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256115.006. 

The  small  version  of  the  problem  is  not  hard  to  solve  even  on  an  ordinary 
computer.  However,  the  large  version  is  more  demanding;  it  required  about  4 
minutes  CPU  time  and  2.7  million  words  computing  storage  space  for  running  one 
single  iteration.  No  more  iterations  were  needed  because  the  solution,  at  the 
end  of  the  first  iteration,  was  the  same  as  the  solution  obtained  by  the  small 
version  of  the  problem  at  that  point;  which  suggests  that  our  choice  for  the 
attractive  routes  was  successful. 

Table  1  gives  a  summary  of  the  intermediate  and  final  results  related  to  the 
solution  of  the  real  world  problem;  the  complete  presentation  of  those  results 
is  provided  by  Appendix  E  in  Ben-Ayed  [1988].  The  Table  includes  only  the 
iterations  and  the  steps  resulting  in  an  improvement  of  the  gap  (increase  of 
the  ideal  solution  or  decrease  of  the  incumbent  solution).  It  can  be  seen  from 
the  table  that,  as  is  usually  the  case  for  iterative  algorithms,  most  progress 
is  made  during  the  first  iterations  while  convergence  is  slower  as  the  number 
of  iterations  increases.  The  first  two  columns  of  the  Table  are  very  different 
from  the  others  because  of  the  instability  of  the  "maximum"  variables  W  in  the 
absence  of  the  corresponding  lower  objective  function;  the  negative  value  of 
the  improvement  in  the  first  column  is  a  consequence  of  this  instability.  The 
fourth  row  of  the  Table  specifies  whether  a  solution  is  ideal  or  incumbent; 
when  a  solution  is  incumbent,  the  rows  5,  6  and  7  refer  to  the  actual  cost 
incurred  by  the  system  (A1),  and  when  the  solution  is  ideal,  those  rows  refer 
just  to  the  value  of  the  objective  function  (01).  Starting  from  the  third 
column,  the  solution  has  a  tendency  to  stabilize;  the  numbers  of  used  links, 
used  routes  and  improved  links  are  almost  the  same  in  all  six  last  columns. 

Several  results  can  be  inferred  from  the  computational  experience  and  the 
solution  of  the  real  world  problem.  Concerning  the  tradeoff  between  efficiency 
and  accuracy,  the  equivalence  of  the  solutions  of  the  large  version  and  the 
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small  version  shows  that  a  complex  problem  can  be  simplified  without  affecting 
the  optimality  of  the  solution.  However  the  simplification  must  be  done 
carefully  and  efficiently;  other  simplifications,  such  as  reducing  the  problem 
to  a  single-level  could  lead  to  bad  results. 

The  actual  cost  incurred  by  the  system  when  the  problem  is  solved  as  an  LP 
exceeds  the  actual  cost  incurred  when  a  BLP  solution  is  used  by  922.224 
thousands  of  Tunisian  Dinars.  There  are  two  reasons  for  not  having  this 
difference  as  high  as  it  is  supposed  to  be.  First,  the  travel  cost  at  stable 
flow  is  a  one-piece  linear  function  that  does  not  depend  on  the  added  capacity. 
However,  if  stable  flow,  which  is  the  overwhelming  proportion  of  the  flow,  does 
not  depend  on  the  added  capacity,  which  is  the  only  instrument  in  the  model  for 
the  system  to  influence  the  user's  choice  of  the  routes,  the  ability  of  the 
system  to  affect  the  user-optimized  equilibrium  is  almost  paralyzed.  Second, 
the  limitation  of  the  added  capacity  resulted  in  low  improvement  costs.  It  can 
be  seen  from  Table  1  that  improvement  costs  are  about  fifty  times  smaller  than 
travel  costs.  The  introduction  of  the  sophisticated  BLP  formulation  of  the 
improvement  costs  did  not  give  much  saving  as  compared  to  the  simpler  LP 
formulation  because  those  costs  are  already  low  according  to  the  data  of  the 
problem.  The  ability  of  BLP  to  give  better  results  is  conditioned  by  the 
availability  of  realistic  data;  providing  more  representative  data  to  our 
problem  would  definitely  lead  to  a  larger  difference  between  the  solution  of 
the  problem  when  solved  as  a  BLP  and  its  solution  when  solved  as  a  single 
level . 

Another  harmful  simplification  would  be  the  inclusion  of  strict  capacities 
on  the  links,  even  when  the  overflow  is  very  small  on  all  links,  as  is  the  case 
in  our  solution.  Imposing  strict  capacities  results  in  the  violation  of 
Wardrop's  first  condition  and  therefore  the  collapse  of  the  equilibrium,  which 
is  one  of  the  bases  of  the  problem. 
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Summary  of  the  Results 

Iteration 

0 

0 

0 

1 

1 

1 

1 

4  &  5 

Step 

- 

- 

- 

1 

2 

3 

4 

4 

Value  of  t 

1.000 

1.000 

1.000 

1.000 

1.000 

1.000 

1.000 

<  .125 

Solution 

ideal 

incumb. 

incumb. 

ideal 

ideal 

incumb. 

incumb . 

incumb. 

Travel  Cost 

255328. 

273736. 

252428. 

245152. 

245320. 

252577. 

251724. 

251680. 

Improvement 

-124325 

3153. 

4609. 

4329. 

4410. 

4410. 

4424. 

4435. 

Total  Cost 

131003. 

276889. 

257037. 

249481. 

249730. 

256987. 

256148. 

256115. 

Routes  Used 

278 

279 

276 

275 

276 

276 

275 

275 

Links  Used 

99 

100 

97 

97 

97 

97 

97 

97 

Linkslmproved 

18 

18 

14 

14 

14 

14 

14 

15 

Gap  (%) 

- 

111.36 

96.21 

3.03 

2.93 

2.91 

2.57 

2.56 

Table  1 


In  the  final  solution,  275  routes  and  97  links  are  used  and  15  links  need 
to  be  improved.  Among  the  112  links  included  in  the  empirical  study,  15  are  not 
used  by  inter-regional  traffic;  one  cannot  tell  whether  those  links  are 
redundant  in  the  inter-regional  network  or  their  lack  of  use  is  due  to  low 
values  of  the  corresponding  entries  in  the  trip  matrix.  The  study  is  more 
concerned  with  the  links  to  be  improved;  those  are  shown  on  the  map  in  Figure 
1.  Almost  all  improvements  happen  to  be  on  the  roads  connecting  Tunis,  the 
capital,  to  its  major  neighboring  cities.  This  result  can  be  intuitively 
predicted  by  a  simple  look  at  the  trip  matrix  [Ben-Ayed,  Blair  and  Boyce  1988]; 
the  trips  originating  or  ending  in  Tunis  are  56  %  of  the  sum  of  all  the  trips 
of  the  matrix.  A  similar  recommendation  was  made  by  the  Plan  Directeur  Routier 
[  SETEC  1982  pp.  6-8]  stating  that  the  roads  corresponding  to  the  exists  of 
Tunis  will  be  highly  congested  by  the  year  2000  if  no  improvements  are  made.  In 
contrast,  the  entries  of  the  trip  matrix  corresponding  to  regions  such  as  the 
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southern  ones  are  very  low,  which  results  in  no  improvement  even  for  unpaved 
roads . 


The 

Links  to  B< 

2  Improved 

Link 

Flow 

Existing 

Added 

Improvement 

Limit  on  Added 

Capacity 

Capacity 

Cost 

Capacity 

1 

2144 

1370 

774 

592 

887 

3 

2537 

1637 

899 

394 

899 

4 

2112 

1115 

997 

297 

1083 

6 

2015 

1055 

960 

615 

1070 

10 

1980 

814 

1166 

484 

1460 

12 

1052 

448 

604 

84 

2046 

14 

306 

271 

35 

22 

1354 

19 

1056 

447 

609 

89 

2046 

20 

2275 

1173 

1101 

463 

1101 

26 

2537 

1637 

899 

245 

899 

28 

2010 

1370 

640 

362 

887 

30 

1517 

1370 

147 

60 

887 

73 

2112 

862 

1250 

394 

1496 

77 

1772 

1055 

717 

156 

1070 

90 

2058 

1370 

687 

179 

887 

Table  2 


In  considering  the  specific  link  improvement,  Table  2  gives  for  each  link  to 
be  improved  the  interregional  flow,  the  existing  capacity  available  to  this 
flow  (60  %  of  the  total  existing  capacity),  the  added  capacity,  its  cost  and 
the  maximum  added  capacity  allowed  by  the  formulation.  The  improvement  required 
for  the  link  connecting  nodes  4  and  29  is  very  small;  a  slight  improvement  of 
the  surface  of  this  link  gives  it  exactly  the  same  capacity  as  the  link 
connecting  nodes  29  and  28  which  does  not  need  improvement  although  it  is 
closer  to  the  capital.  All  other  improvements  are  more  significant  and  three  of 
them,  namely  those  on  links  3,  20  and  26,  require  the  maximum  added  capacity 
allowed  by  the  formulation.  In  each  of  those  links  the  flow  is  higher  than  the 
new  capacity;  more  detailed  study  is  needed  to  include  the  possibility  of 
upgrading  them  to  three  or  four  lane  highways.  Link  90  is  also  congested; 
however  the  decision  maker  responsible  for  the  system  does  not  add  more 
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capacity  although  he  has  the  possibility  to  do  so.  The  reason  is  that  the 
congestion  on  those  roads  is  less  expensive  than  the  addition  of  more  capacity. 

4.  Conclusions 

BLP  remains  a  little  known  technique,  rarely  used  despite  all  its  potential 
capabilities.  The  computational  difficulties  of  BLP  may  not  be  the  main  reason 
of  its  infrequent  applications;  there  are  other  optimization  techniques,  such 
as  Mixed  Integer  Programming,  that  are  not  easier  to  solve  but  are  widely  used. 
BLP  is  a  new  model;  a  great  effort  is  needed  to  help  people  know  it  and  take 
advantage  of  its  powerful  formulations.  Developing  efficient  and  reliable 
software  that  uses  data  structure  and  that  includes  several  heuristics  and 
branch  and  bound  algorithms  would  have  considerable  effect  on  increasing  the 
use  of  BLP. 

The  empirical  study  conducted  in  this  research  demonstrated  that  BLP  is 
tractable  and  can  be  efficiently  used  to  formulate  and  solve  even  complex  real 
world  problems  such  as  the  design  and  the  improvement  of  a  transportation 
network.  The  complexity  of  BLP  should  not  be  a  reason  for  formulating  bilevel 
problems  as  single- levels ;  the  solution  to  the  single-level  problem  might 
happen  to  be  far  away  from  the  actual  optimum,  and  the  gain  in  simplicity  and 
computation  time  may  not  be  worth  the  loss  in  optimality.  By  taking  into 
account  the  special  structure  of  the  problem  at  hand  and  making  use  of  the 
powerful  computers  available  to  us  nowadays,  it  is  always  possible  to  solve  the 
BLP  problem  and  get  a  better  solution  than  the  one  obtained  when  approximating 
it  by  a  single  level.  For  any  BLP  problem,  the  use  of  any  related  algorithm 
always  gives  a  solution  that  is  better  than,  or  at  worst  the  same  as,  the 
solution  obtained  when  ignoring  the  lower  objective  and  solving  the  problem  as 
a  single-level.  In  other  words,  the  relative  inefficiency  of  BLP  algorithms 
does  not  justify  the  formulation  of  bilevel  problems  as  single-levels. 
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Appendix 

A.    Definition  of   the  Variables    and  the  Functions    in  Formulation   (2) 

F(Z,W)    =   Ea=liU2   Fa(Za>Wma)   =   Ia=ljll2    ((blaZa+b0a)    +   Em=1>Ma(bm+lja-bma)Wma 

F(W)   =   ^a=ljii2  Ia(Wma)  =   Za=l,112    Zm=l,Ma  Wma 


W)   =  Za=l,112   GaCCa) 
G(C)    =   Ea=l,U2   Ga(Ca) 

H1(Y)  =  0  is  defined  as: 
H2(X,Y)  =  0  is  defined  as: 
H3(Z,C,X)  >  0  is  defined  as: 


=  Za=l,112  ^Ca--4c0aka) 
=  Za=l,112  ^a 


H4(Z,C,X)  >  0  is  defined  as: 


H5(Z,W)  >  0  is  defined  as: 


H^(Z)  >  0  is  defined  as: 


for  each  origin-destination  pair  od: 


ZreRod  Yr  ~  uod 
for  each  link  a: 


Xa  "  Zr=l,2095  5arYr  ~  ° 
for  each  link  a: 

Ca  "  <WXa  +  -4ka)  *  ° 

ca  "  cla(xa  +  -4ka)  +  (cla  "  <Wza  +  dla  *  ° 
for  each  link  a: 

Ca  *  c0a<xa  +  -4ka) 

Ca  *  ^la^a  +    -4ka)    "    tela    "   ^Oa)Za  +  <ha 
for   each    link   a: 

Wma   "    Za  +  %,a  ^   °  m=1>Ma 

for   each    link   a: 


Z„  ■-  q3r  >  0 


'a   ^3a 
Z,  C,  C,  W,  X,  Y  >  0  is  defined  as:  for  each  link  a: 

Za>  Xa'  Ca'  ^a>  Wma  -  ° 


m=l,M. 


for  all  reR 


od: 


Yr  >  0 


where  Cq3  and  c^a  are  the  slopes  of  the  system  travel  cost  on  link  a  at  stable 
and  unstable  flows,  respectively;  cQa  and  cla  are  the  slopes  of  the  cumulative 
user  travel  cost  at  stable  and  unstable  flows,  respectively;  dla  and  dla  are 
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the  intercepts  of  the  travel  cost  at  unstable  flow  of  the  system  and  the 
cumulative  user,  respectively;  Mfl  is  the  number  of  breakpoints  in  the  improve- 
ment cost  function  in  link  a;  Mfl  is  0  when  the  curve  is  linear,  1  when  the 
curve  has  2  pieces  and  2  when  it  has  3  pieces;  bga  is  the  intercept  of  the 
first  piece  of  the  improvement  cost  function  on  link  a;  bma  is  the  slope  of  the 
piece  delimited  by  qm_^  a  and  qma,  m=l,Ma+l,  qQa  being  equal  to  zero  and  q^a+^ 
equal  to  q3a,  the  maximum  improvement  allowed  by  the  model;  kfl  is  the  existing 
capacity  of  link  a;  Xa  is  the  flow  on  link  a;  Cfl  is  the  system  total  travel 
cost  on  link  a;  C_a  is  the  cumulative  user  total  travel  cost  on  link  a;  Za  is 
the  number  of  passenger  car  units  (PCU)  added  to  the  capacity  of  link  a;  W_a  is 
the  maximum  of  (Za~qma)  and  0;  Y  is  the  flow  on  route  r;  R  j  is  the  set  of  all 
routes  from  origin  o  to  destination  d;  uQ(^  are  the  entries  of  the  trip  matrix; 
6ar  is  binary  number  equal  to  1  when  link  a  belongs  to  route  r,  and  equal  to  0 
otherwise. 

B.  Description  of  the  Algorithm 

When  describing  the  algorithm,  we  use  T,  T   and  T  ,  01  as  defined  in  Section  2: 

Iteration  0 

•  Find  a  current  ideal  solution,  a  current  incumbent  solution  and  calculate 
the  current  gap  for  BLP  problem  (2). 

•  Call  Termination  test. 

•  Set:   i  =  1,   xMAX  =  1,   tmin  =0,   1=1. 

•  Start  iteration  I. 
Iteration  I 

Set:   T1  =  <f>,      T2  =  T. 
Step  1 

•  Solve  the  LP  obtained  after  substituting  all  nonlinear  improvement 
functions  by  aZ  +  6  for  all  links  to  find  the  links  of  T  that  are 
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candidates  for  improvement. 

•  Get  the  solution  (Z1,  X1,  C1,  C1,  Y1,  W1). 

•  For  every  link  j  of  T  ,  if  the  corresponding  Zi       is  not  equal  to  zero, 

7  1 

subtract  j  from  T  and  add  it  to  T  . 

Step  2 

•  Solve  the  BLP  problem  (7)  with  the  set  of  lower  variables  equal  to  {Wma: 

bm+l,a-bma  <  °>  m=1>25  asTl)' 

•  Get  the  solution  (Z2,  X2,  C2,  C2,  Y2,  W2). 

•  If  there  are  elements  j  of  T  such  that  Z-   is  positive,  then: 

9  i 

Subtract  from  T   all  elements  j  and  add  them  to  Tx 
Repeat  Step  2 
Else: 

The  optimal  solution  of  (6)  is  obtained 

End  If. 

2 

•  If  x  =  1  and  0   >  value  of  current  ideal  solution,  then: 

Update  ideal  solution  and  calculate  new  gap 
Call  termination  test 
End  If. 
Step  3 

•  Solve  the  lower  LP  problem  (5)  with  Z  equal  to  Z. 

•  Get  the  solution  (X3,  C3,  C3,  Y3) . 

9 

•  If  A   <  value  of  current  incumbent  solution,  then: 

Update  incumbent  solutions  and  calculate  new  gap 

Call  termination  test 
End  If. 

•  If  A2  =  02  (or  A2  ~  02),  then: 

(Z2,  X2,  C2,  C2,  Y2,  W2)  is  feasible 


Set:   x 


MIN  "  T>   x  =  (TMAX  +  T)/2 
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Else: 

Try  to  improve  feasible  solution 
Set:   xMAX  =  T,   t  =  (xMIN  +  t)/2 
End  If. 
•1  =  1  +  1. 
•  Start  iteration  I. 
Termination  Test 

•  If  I  >  N,  then: 

Number  of  iterations  exceeded  before  reaching  desired  accuracy 

Incumbent  is  the  best  solution  found 

Stop 
Else  if  gap  <  e: 

Incumbent  solution  is  optimal  within  the  specified  accuracy  e 

Stop 
End  If. 

•  Return. 

Try  to  improve  feasible  solution 

•  Solve  (6)  while  locking  the  nonbasic  variables  (including  slack  and 
surplus  variables)  of  the  LP  problem  (5)  solved  with  Z  =  Z  . 

•  Get  the  solution  (Z4,  X4,  C4,  C4,  Y4 ,  W4). 

•  If  0   <  value  of  current  incumbent  solution,  then: 

Update  incumbent  solution  and  calculate  new  gap 
Call  termination  test 
End  If. 
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